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Abstract 

This  paper  addresses  the  formulation  of  hierarchical  plate  models  as  an  op¬ 
timal  (in  a  clearly  defined  sense)  numerical  method  which  can  be  easily  imple¬ 
mented  and  is  available  in  the  code  MSC/PROBE.  Special  emphasis  is  placed  on 
the  optimal  selection  of  the  shear  correction  factors.  It  is  shown  that  different 
measures  of  accuracy  lead  to  different  optimal  choices  of  these  factors.  The  main 
tool  in  the  analysis  is  the  Fourier  transformation. 


1.  Introduction 

The  problem  of  plate  modelling  is  a  classical  problem  of  mechanics.  The  principles 
underlying  the  derivation  of  plate  models  can  be  divided  into  .3  broad  groups. 

a)  Physical  derivation.  Here  assumptions  of  a  geometrical  nature  are  made  a  priori 
in  conjunction  with  additional  assumptions  involving  certain  stress  components. 
A  typical  example  is  the  KirchhofF  hypothesis  [26]  which  leads  to  the  well  known 
biharmonic  equation  for  the  normal  deflection  of  the  plate.  This  approach  has 
been  generalized  in  many  papers  and  we  refer  here  especially  to  [1,  2,  10,  17,  23, 
24,  25,  30,  31,  34,  35,  37,  39,  40,  41,  44,  46,  49]  and  the  references  therein.  The 
idea  of  a  shear  correction  factor  was  first  introduced  in  the  papers  of  this  group. 

b)  Asymptotic  Analysis.  Here  plate  models  are  justified  as  the  leading  terms  of  an 
asymptotic  expansion  of  the  three-dimensional  solution  with  respect  to  the  plate- 
thickness.  We  refer  here  for  example  to  [19],  [20]  and  to  [34],pp.  585-588,  w'here  a 
historical  account  for  this  approach  is  given.  Here  the  asymptotic  methods  were 
mostly  applied  to  the  plate  problem  formulated  as  a  differential  equation  rather 
than  to  its  variational  form.  Asymptotic  methods,  where  a  mixed  variational 
principle  was  used,  appeared  only  recently  [11].  See  also  [7],  [33]  for  related  ideas. 

c)  Numerical  and  Hierarchical  Approach.  Here  we  understand  the  modelling  as  ap¬ 
proximate  solution  of  the  three  dimensional  problem  (which  is  solved  for  example 
by  the  finite  element  method).  It  is  essential  to  use  an  approach  which  can  be 
easily  implemented  (  and  which  avoids  in  particular  higher  order  equations)  in  the 
framework  of  the  h-p  version  of  the  finite  element  method  and  which  has  optimal 
properties  yielding  the  most  accurate  results  when  compared  with  the  exact  solu¬ 
tion  of  the  three  dimensional  problem.  It  has  to  be  emphasized  that  various  plate 
models  can  significantly  differ  in  the  area  where  boundary  layers  and  singularities 
dominate  the  solution  behaviour.  Hence  a  flexible  hierarchy  is  essential  (and  has 
been  implemented  in  the  code  MSC/PROBE).  A  posteriori  error  estimators  and 
adaptive  approaches  for  the  model  selection  within  the  hierarchy  are  yet  to  be 
developed. 


In  this  paper  we  concentrate  on  the  selection  of  the  shear  correction  factor  which 
leads  to  the  highest  accuracy  of  the  solution.  We  restrict  ourselves  to  the  case  of  an 
infinite  plate  made  of  a  homogeneous,  isotropic  material  and  use  the  Fourier  transform 
technique  in  our  analysis.  The  problem  of  anisotropic  and  laminated  plates  will  be 
discussed  elsewhere.  In  the  second  section  we  define  the  three  dimensional  problem 
and  the  hierarchy  of  the  models.  The  third  section  addresses  the  Fourier  transform 
of  the  solution  for  the  three  dimensional  formulation  as  well  as  for  the  plate  model 
problem.  In  section  four  we  analyze  the  errors  of  various  models  and  give  the  optimal 
shear  fa^ctors.  In  section  five  we  present  numerical  examples  and  in  section  six  we  sum 
up  our  basic  conclusions. 


2.  Hierarchical  Plate  Models 

2.1.  The  Plate  Problem 

The  plate  problem  is  the  boundary  value  problem  of  three  dimensional,  linear  elasticity 
for  an  isotropic  material  on  the  domain 


=  {x  =  (xi,X2,X3)  I  (xi,X2)  €  w,  ksl  <  d}  =a;  X  {-d,d)  (2.1) 


where  w  C  is  a  bounded  domain  with  a  piecewise  smooth  boundary.  Here  u)  is 
referred  to  as  the  midsurface  and  2d  as  the  thickness  of  the  plate. 

As  usually,  we  denote  by  u  =  {u,},  i  =  1,2,3  the  displacement  vector.  By  cr  = 
Wij}i  hj  =  1,2,3  and  e  =  {e,j},  i,j  =  1,2,3  we  denote  the  stress  and  strain  tensor, 


respectively,  with  e 
and  strain  can  then 


be  written  in  the  form 


1,2,3.  Hooke’s  law  relating  the  stress 
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^  “  2(1  +J/)  -  <^  =  =  «5S  =  “66, 

and  a,j  =  0,4  <  i  <  6,  j  =  1,...6,  j  ^  i. 

Here  E  is  Young’s  modulus  and  v  is  the  Poisson  ratio. 

We  will  focus  on  the  classical  plate  bending  problem  where  the  normal  tractions  ate 
prescribed  on  the  upper  and  lower  plate  surface.  For  (ii,X2)  €  a;  we  have 

<^33(2:1,  X2,  ±d)  =  ±i9(Xi,X2) 

<7’a3(xi,X2,±d)  =  0,  a  =  1,2 

where  9(xi,X2)  €  L2{u}).  The  solution  u  is  the  minimizer  of  the  total  plate  ener<iv 
Gp{u)  =  -f^(u)  -  ^ -g(x)(u3(xi,X2,<f) +  U3(xi,X2,-d))dxidx2,  (2.6  I 

where 

S^iu)  =  /  e-^Aedx  =  ||u||| 

JQP 

is  the  strain  energy,  over  a  suitable  space  'H(Q'^)  C  where  characterizes 

the  boundary  conditions  on  the  lateral  plate  sides  and  H^{Vl)  denotes  the  usual  Sobolex 
space. 

By 

=  {u  I  u  =  (01,02,03)  +  {bi,b2,b^)  X  (xi,i2,X3),  ai,bi  ^  R  i  =  1,2,3}  (2.7) 

we  denote  the  set  of  rigid  body  motions.  Now,  if  A^(np)  =  n  =  0  and 

q  G  ^2(0;),  the  quadratic  functional  Qp{u)  has  a  unique  minimizer.  This  follows  from 
Korn’s  inequality  (for  an  elementary  proof  of  which  we  refer  to  [27])  and  the  standard 
use  of  the  Lax-Milgram  lemma. 

It  is  obvious  that  due  to  the  symmetry  of  the  problem,  we  can  constrain  the  space 
of  solutions  to  be  of  the  form 

Ui(Xi,X2,X3)  =  _Ui(li,X2,-X3) 

«2(a;i,X2,X3)  =  -«2(Xl,X2,-X3)  (2.S) 

«3(ll,X2,a:3)  =  U3(Xi,X2,  -X3). 

Relation  (2.8)  constraints  the  space  accordingly  and  will  be  assumed  everywhere 

in  what  follows.  In  the  case  that  H  K(n'’)  =  ^  0  we  have  to  assume  in 

addition  that 

Jj0,0,q).n^dx  =  0  V7e€W^(np),  (2.;)) 
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then  the  minimizer  of  (2.6)  exists  and  is  unique  modulo 

For  w  =  and  :=  u  x  (— d,d)  the  situation  is  slightly  more  complicated.  In 
this  case  we  will  assume  that  g(xi,X2)  €  L2{R^),  has  compact  support  and  satisfies 
(2.9)  for  any  1Z  €  N’’’.  Then  the  quadratic  functional  Qp  has  a  minimizer  over 
where 

n^lo)  =  K  €  i  =  1,2,3  I  f  e^Aedx  <  oo}.  (2.10) 

•/nSo 

The  space  'H(n^)  is  not  a  subspace  of  because  the  displacement  components 

Ui  €  'H(n^)  are  in  general  not  square  integrable.  For  more  about  the  analysis  of  the 
characteristics  of  spaces  of  the  type  'W(n^)  and  related  questions  we  refer  to  [27].  The 
minimizer  of  is  unique  modulo  functions  in  and  once  more  we  will  constrain  the 
space  'H(fl^)  by  (2.8). 

We  will  also  address  the  q  =  (ai, a!2)-periodic  plate  problem  on 

fip  =  Wa  X  {-d,d),  Wa  =  (0  <  Xi  <  Q,-,i  =  1,2). 

Here  we  will  consider  spaces  of  a-periodic  functions  in  (xi,  X2)  with  periods  in  xi,  .i  ;. 
A  function  <^(xi,  X2,  X3)  is  called  a-periodic  if  <i>{xi  +  ma\,X2  +  na2,  X3)  =  <?!>(xi,  X2,  .f3 ) 
for  arbitrary  integers  m,  n.  We  define 

-W{a)  =  {u,  €  (.Him),  i  =  1,2,3)  (2.11) 

where  is  the  space  of  a-periodic  functions  in  (xi,X2)  and  analogously  we  also 

define  L2,a{^a)-  The  minimizer  of  Q^^q{u)  is  now  defined  by  (2.6)  with  the  integration 
taken  over  and  Wa,  respectively.  The  solution  of  the  o-periodic  plate  problem  is  the 
minimizer  of  ^^a(^)  ^^(^)-  This  minimizer  exists  if  the  q  satisfies  (2.9)  (where 

the  integral  is  now  taken  over  Wa)  for  all  Tt  €  fl  'H{a)  =  and  is  unique  modulo 
AfP.  Once  more  the  constraint  (2.8)  will  be  assumed  here. 

Together  with  the  plate  problem  we  will  also  consider  its  2-dimensional  counterpart 
which  we  will  refer  to  as  the  beam  problem.  Its  formulation  can  be  obtained  formall\ 
from  that  of  the  plate  problem  by  the  assumptions 


In  this  case  we  have 


9(xi,X2)  =  q{xi) 


dui 

5X2 


=  U2  =  0,  i  =  1,3. 


n**  =  {x  =  (x,,X3)  I  |xi|  <  L,|x3|  <  d}  =  w  X  (-d,d) 


(2.12) 


0 


where  <*;  =  /  =  {—L,L).  The  displacement  vector  is  now  u  =  («i,U3),  <t  =  Wij], 
e  =  {cij},  i,j  =  1,3  and  Hooke’s  law  is  given  by 

<Tn  A  +  2/i  A  0  cii 

(T22  —  A  A  +  2/i  0  =  633  (2.13) 

<T33  0  0  /i  2ei3 

or 

<r  =  Ae.  (2.14) 

As  before,  for  q  ^  L2(I),  u  is  the  minimizer  of  the  quadratic  functional 

^  £qixi)(u3{xi,d)  +  usixi,  -d))dxi  (2.15) 

where  S^{u)  =  f^b  e^Aedx  =  ||u||£  over  a  suitable  space  C  satisfying 

the  constraint  (2.8)  for  Ui,  U3  which  characterizes  the  boundary  conditions  on  the  lateral 
sides.  We  denote  as  before  by 

W*’  =  {^(11,12)  I  =  (<11,03)  +  6(x3,-Xi),  ai,be  Ri  =  1,3}  (2.16) 

the  set  of  rigid  body  motions.  Then  if  A/’^fob)  =  A/’*’  H  =  0,  the  solution  exists 

and  is  unique.  Otherwise  we  have  to  assume  that 

J^{0,q)-n^dxi  =  0  (2.17) 

In  this  case  the  minimizer  exists  and  is  unique  modulo  Af^^^by 

U  (jj  =  and  17^  =  R^  x  {—d,d)  then  we  assume  that  ^(xi)  G  L2{R)  hzis  compact 
support  and  satisfies  (2.17)  for  any  TZ.  £  Af^.  The  minimizer  of  over 

W(Hto)  =  {o.  €  {HU^D),  i  =  1,3  I  f  e'^Aedx  <  00}  (2.1S) 

exists  and  is  unique  modulo  Af’  as  for  the  plate  problem.  We  note  that  7f(n^)  is  not 
a  subspace  of  (/f'(H*’))^  because  u,-,!  =  1,3  do  not  have  to  be  square  integrable  on  9A. 
Quite  analogously  as  above  we  define  the  a-periodic  beam  problem.  Here  we  set 

=  Wo  X  {—d,d).  Wo  =  (0  <  ii  <  a) 

and  will  consider  the  space  of  a  periodic  functions  (in  ii)  on  Ho-  define  as  before 

7f''(a)  =  {u.€//^(n‘),f  =  l,3}  (2.19) 

where  space  of  a-periodic  functions  in  xi  with  the  condition  (2.8).  Analo¬ 

gously  we  define  L2,flr(<^o).  The  functional  ^^o(^)  defined  by  (2.15)  where  the  integrals 
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are  now  taken  over  D*  and  Wq,  respectively.  Then  the  solution  of  the  a-periodic  beam 
problem  is  the  minimizer  of  This  minimizer  exists  if  q  satisfies 

(2.17)  (where  the  integral  is  to  be  taken  over  Wq)  for  all  71  €  M’’’’  H  and  it 

is  unique  modulo 

Here  and  in  what  follows,  we  are  using  the  super-  and  subscripts  p  and  h  to  underline 
that  we  have  in  mind  the  plate  or  beam  problem,  respectively.  If  no  ambiguity  could 
occur,  however,  we  shall  simply  omit  these  super-  and  subscripts. 

2.2.  Hierarchical  Modelling  of  Plates  and  Beams 

The  fact  that  a  parameter  d  which  is  small  with  respect  to  the  diameter  of  u  suggests 
to  replace  the  three-dimensional  plate  problem  by  a  simpler,  two-dimensional  plate 
model.  Since  every  plate  model  is  only  an  approximation  of  the  three-dimensional  plate 
problem  our  goal  is  to  construct  a  hierarchy  of  plate  models  (of  increasing  complexity) 
that  are  capable  of  approximating  the  solution  of  the  three-dimensional  plate  problem 
for  any  fixed  thickness  d  arbitrarily  closely. 

To  this  end  we  first  introduce  the  notion  of  the  n-model. 

Definition  2.1  Let  n  =  (n,),i  =  1,2,3  he  a  vector  of  nonnegative  integers  and  B{n)  = 
{6ij}  be  a  (Hooke’s)  matrix  of  the  type  (2.2)  which  is  assumed  to  be  symmetric  and 
positive  definite  and  in  general  5(n)  ^  A.  Then  we  denote  by  ^“(x)  the  solution  of  the 
n-plate  model  which  is  the  unique  minimizer  ofQp{u)  over  the  subspace  7i{n)  C 
of  displacement  fields  of  the  form 


=  ■  (2.20) 

j=o  \  a  / 

For  the  beam  the  definition  is  analogous. 

Existence  and  uniqueness  of  the  solution  of  the  n-model  are  proved  along  the  same 
lines  as  for  the  exact  solution.  We  can  once  more  consider  the  bounded  domain  u,-  or 
the  entire  space  H}  for  the  plate  (respectively  R}  for  the  beam). 

Remark  2.1  The  (material)  matrix  B  used  in  the  definition  of  the  n-model  may,  and, 
in  fact  must  be  at  times  different  from  the  matrix  A  of  the  plate  problem,  as  we  shall 
show  below. 

Remark  2.2  The  approximation  (2.20)  is  a  classical  one.  The  case  n  =  (1, 1,0)  was 
used  for  example  in  [10],  [14],  [40],  n  =  (1,1,2)  in  [17],  [49],  n  =  (2,2,2)  in  [35]  and 
n  =  (3, 3, 2)  in  [31,  41].  In  addition  these  forms  have  been  used  when  a  relation  between 
Uij  (by  differentiation)  is  imposed  a  priori,  as  e.g.  in  [46].  Many  different  approaches  to 
derive  the  differential  equations  describing  the  model  were  suggested.  We  refer  here  foi 
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example  to  [1,  2,  23,  24,  25,  30,  34,  37,  44],  The  form  (2.20)  was  also  used  in  conjunction 
with  various  variational  principles  (direct  or  complementary)  and  additional  heuristic 
ideas  to  obtain  the  differential  equations  describing  the  model.  The  exact  solution  to 
the  three-dimensional  plate  problem  and  the  error  estimation  have  not  been  considered 
in  these  papers. 

For  the  asymptotic  error  analysis  called  “justification”  of  various  models  we  refer 
to  [11,  12,  14,  15]  and  for  an  assessment  and  classification  of  various  approaches  as 
well  as  for  additional  references,  see  for  example  [12],  and  [7,  8,  33,  34].  The  main 
question  here  is  the  analysis  of  the  difference  between  the  solution  of  the  model  and 
three  dimensional  solution  as  d  —*  0.  In  these  papers  mainly  the  two  sided  estimates 
based  on  direct  and  complementary  variational  principles  are  used. 

In  the  papers  [9,  10,  16,  23,  24,  25,  30,  35,  36,  39,  49,  50]  the  form  (2.20)  is  also 
suggested  for  laminated  plates.  In  this  context  the  following  question  arises: 

Is  the  form  (2.20)  a  heuristic  one  or  can  it  be  rigorously  justified? 

In  the  case  of  homogenous,  isotropic  material  it  wcis  shown  in  [45]  that  (2.20)  is  not 
optimal  for  d  >  0  and  that  a  smaller  error  with  respect  to  the  three  dimensional 
solution  can  be  obtained  with  a  model  of  the  same  complexity.  Nevertheless,  (2.20)  is 
asymptotically  optimal  as  d  — +  0.  For  laminated  plates,  however,  (2.20)  is  neither  for 
d  >  0  nor  in  any  asymptotic  sense  optimal.  An  optimal  form  exists  when  the  shape 
functions  in  X3  belong  to  the  kernel  of  a  differential  operator  intrinsic  to  the  problem, 
as  hcis  been  shown  in  [47],  [48].  For  a  concrete  example  of  the  hierarchical  modelling 
of  laminated  plates  along  these  lines  we  refer  to  [6]. 

Our  approach  uses  the  minimization  of  the  total  energy,  i.e.  the  approximate  solu¬ 
tion  is  a  projection.  The  merit  of  this  method  lies  in  that  it  can  be  very  naturally  and 
efficiently  implemented  within  the  frame  work  of  the  h-p  version  of  the  finite  element 
method  which  is  available  in  the  code  MSC/PROBE.  Hence  the  hierarchical  modelling 
has  to  be  understood  as  a  numerical  method  for  the  solution  of  the  original  3  dimen¬ 
sional  problem  with  its  singularities  and  boundary  layers.  The  a-priori  and  a-posteriori 
error  estimation  eis  well  as  the  adaptive  model  selection  are  essential  here. 

Since  we  consider  only  the  problem  of  pure  bending,  the  exact  solution  of  the 
plate  problem  satisfies  (2.18).  Correspondingly,  we  may  set  in  (2.20) 

Ua,j  =  {x\,X2)  =  ^  for  ;■  even  ,Q  =  1,2,  . 

f/3j(xi,i2)  =  0  for  j  odd. 

Since  ^“(x)  is  only  an  approximation  to  the  solution  uP(x)  of  the  plate  problem, 
we  will  judge  the  accuracy  of  the  n-model  by  the  relative  error  in  a  quantity  of  interest 
when  a  load  q  \s  given.  Here  'P  is  a  suitable  class  of  loads  which  will  be  specified 
later. 
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For  example,  if  we  are  interested  in  the  error  in  energy,  we  define 


E{d,n,q) 


I  5^  (up)  ) 


(2.22) 


where  S‘^(u’’)  resp  is  the  strain  energy  of  u’’  respectively  of  Here  V 

could  be  the  class  of  loads  with  compact  support  which  satisfy  (2.9).  Then  (2.22)  is 
obviously  well  defined  for  any  q  €.V. 

Alternatively  to  (2.22),  we  can  also  consider  the  measure 


E[d,x\,q) 


/|£"‘(^P-“)-5^(uP)|\ 

V  J 


1/2 


(2.23) 


where  could  be  a  “postprocessed”  solution  constructed  from 
For  example  we  may  be  interested  in 


E(d,  n,  q)  = 


(/fli  (^3(2;i,X2,0)  -®u'’’“(ii,X2,0))^dxidx2) 
(fjt3(n§(xi,X2, 0))^dxidx2y^^ 


1/2 


(2.24) 


where  the  class  of  loads  P  is  such  that  for  a  proper  selection  of  the  rigid  body  motion 
the  numerator  and  the  denominator  are  well  defined  for  every  q  £P.  We  remark  that 
the  error  measure  in  (2.22)  -  (2.24)  are  relative  errors.  They  are  natural  in  that  in  their 
analysis  no  assumptions  on  a  scaling  of  the  data  q  has  to  be  made  a  priori.  We  will 
show  below,  however,  that  the  accuracy  of  the  model  can  be  measured  by  analyzing 
only  the  numerator  in  (2.22)  -  (2.24). 

Let  us  now  define  precisely  what  we  mean  by  hierarchical  models  and  give  the 
framework  in  which  we  analyze  the  above  error  measures. 


Definition  2.2  A  sequence  of  n-models  is  hierarchical  with  respect  to  the  measure  E 
and  the  class  P  of  loads  q  if 

(HI)  (asymptotic  exactness)  for  all  q  €P,  E{d,n,q)  —*  0  as  d  0  for  all  n. 

(H2)  ifP  is  such  that  the  exact  solution  is  sufficiently  smooth  uniformly  with  respect 
to  d  then 

E{d,n,q)  <  C{n)d^^^\  /3{Tn)>P{n)  for  >  n^,  z  =  1,2,3. 

The  exponent  ^(n)  will  be  called  rate  of  convergence. 

(H3)  For  any  q  GP  and  any  fixed  d  >  0  there  holds  E{d,n,q)  —*  0  as  min{ni)  — ►  cc, 
1  <  z  <  3. 
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(H4)  The  plate  model  is  mo-iotonically  hierarchical  if  for  sufficiently  largt  n  and  n,  > 
Hi,  i  =  1,2,3,  E{d,m,q)  <  E{d,n,q). 

For  the  beam  problem  the  definition  is  analogous. 

Remark  2.3  It  is  clear  that  in  our  approach  not  only  the  solution  of  the  plate  models, 
but  also  the  modelling  itself  (i.e.  the  selection  of  a  particular  model  from  the  hierarchy) 
could  be  done  numerically.  This  is  in  contrast  to  the  approaches  proposed  in  the 
references  in  Remark  2.2. 

The  following  Theorem  establishes  the  important  properties  (H3)  and  (H4)  for  our 
hierarchical  models. 

Theorem  2.1  Let  the  n-model  be  as  in  Definition  2.1  with  B{n)  =  A  for  all  n.  Then 
for  the  measures  (2.22)  or  (2.23)  and  for  the  conditions  (H3)  and  (H4) 

hold. 

Proof:  is  merely  the  energy  projection  of  onto  'W(n).  This  in.plies  (H4)  lor 

(2.21).  The  density  of  the  polynomials  in  //'*'(— d,d).  A:  =  1,2  yields  (H3)  for  (2.22). 
Realising  that  for  we  have  |£'^(4i’’'“)  —  5"^(u’’)|  =  —  u^),  and  (H3) 

and  (H4)  for  the  measure  (2.23)  follow.  □ 

Remark  2.4  It  is  essential  here  that  we  formulate  our  models  as  energy  projections 
and  not  in  the  form  of  asymptotic  expansions  which  cecise  to  make  any  sense  if  the 
three  dimensional  solution  is  not  smooth.  This,  however,  happens  frequently  because 
of  the  corners  and  edges  that  are  typically  present  and  because  of  the  boundary  layers 
in  the  three  dimensional  solution.  In  this  situation  many  models  eis  e.g.  some  of  those 
in  [2]  are  not  usable. 

We  emphasize  here  that  the  hierarchical  approach  allows  to  employ  simultaneous!}- 
models  of  different  order  n  from  the  hierarchy  in  different  parts  of  the  plate.  In  par¬ 
ticular  in  a  neighborhood  of  the  boundary  of  uj  higher  order  models  (in  our  sense)  are 
essential,  because  among  other  things  boundary  layers  are  different  for  different  models, 
see  e.g.  [4],  [5]. 

Theorem  2.i  insures  the  desirable  properties  (H3)  and  (H4)  for  the  n-models  and 
the  measures  (2.22),  (2.23).  Obviously  these  properties  hold  only  for  large  n  and  the 
proof  sheds  no  light  on  the  properties  of  the  n-models  for  small  n  and  d.  In  practice, 
however,  we  are  especially  interested  in  these  lowest  members  of  the  hierarchy.  These 
are  in  view  of  (2.20),  (2.21)  the  (1,1,0)-,  (1,1,2)-  and  the  (3, 3, 2)-model,  respectively. 
In  the  following  sections  we  therefore  analyze  these  models  in  detail  and  show  that  they 
can  be  optimized  in  some  sense. 
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3.  Fourier  Transform  Solution  of  the  Plate  and  Beam  Problem 
and  Their  n-Models 

3.1.  The  Fourier  Transform  of  the  Plate  and  Beam  Problem 

In  Section  2  we  defined  the  plate  problem  for  u  =  and  the  beam  problem  for  u  = 
as  well  as  the  corresponding  o-periodic  problems. 

In  this  section  we  will  further  elaborate  on  these  two  problems  using  the  Fourier 
transformation  as  the  main  tool.  Because  we  have  to  deal  with  functions  which  are 
not  integrable  we  will  employ  the  theory  of  Fourier  transform  of  generalized  functions 
(tempered  distributions),  see  for  example  [18]. 

For  any  integrable  function  q{x)^  i  we  denote  by 

qiO=  [  9(x)e*^dx  ^  G  i?"  (3.1) 

JR’' 

the  Fourier  transform  of  q{x).  Then  q(x)  is  defined  by  the  inverse  formula 

?(x)  =  (2t)-"  /  (3.2) 

JR” 

If  q{x)  €  L2{R^)  then  9(^)  €  L^iR'^)  and  Parseval’s  formula 


f  (f>^  dx  —  (2x)  ”  ^  p 

JR”  JR” 


holds. 


Tr  \  •  (2-k  2t:  2-k\  .  , 

11  q[x)  is  an  Q  =  I  — ,  — ,  •  •  • ,  -5-  j  periodic  function  and 

\  Pi  P2  Pn ) 


=  (3.4) 

m 

where  [m/3]  =  {mif3i,m2/32,  ■ .  ■  ,mn0n),  rn  =  (mi,  m2, . . . ,  m„),  m,  an  integer,  then  in 
the  sense  of  generalized  functions 


m 

where  <5(^)  is  the  Dirac  distribution. 

The  inverse  transformation  is  then  given  by  formula  (3.4)  and  for  any  r(x)  = 
(27r)~"  Ylm  Parseval’s  formula  reads 

/  gfdx  =  (2T)-"(A/?2.../?n)-‘X^-^-  (3-C) 

JO<Xi<2n/0,  „ 
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If  q{x)  €  L2{R^)  and  has  compact  support  then  q{^)  has  all  derivatives  in  a  neigh¬ 
bourhood  of  ^  =  0  and  the  condition  (2.9)  reads 

«0)  =  ^(0)  =  0  <=I,2.  (3.7) 

For  the  periodic  problem  condition  (2.9)  yields 

9(0)  =  Co  =  0.  (3.8) 


Analogous  results  hold  for  the  beam-problem. 

In  what  follows  we  denote  by  Ui(^i,^2»3:3),  i  =  1,2,3,  the  partial  Fourier  transforms 
of  the  displacement  components  with  respect  to  (xi,X2)  ^^nd  in  the  case  of  a  beam  we 
will  write  analogously  u,'(<fi, ^^3),  with  i  =  1,3. 

Now  using  Parseval’s  formula  we  have 

Theorem  3.1  Let  u^(x)  resp  u*’(x)  have  finite  strain  energy.  Then  for  the  plate 


=  /  e^  Aedx 


2(1  -f- 


-I-  +  Uapj  diid^idxj. 


Here  =  ^1  +  ^2  prime  indicates  - — .  Further,  for  the  beam 

0x3 

+  luil^  +  |u3p  +  ^  _  2j,  d"  I  d^idx^. 
For  the  periodic  case  when  writing 

we  get  for  the  plate 

f^(tif)  =  f  e^  Aedx 

+  -f-  dx^ 


(3.9) 


(3.10) 

(3.11) 

(3.12) 
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where  m  =  (mi, m2).  For  the  beam  we  get  analogously 

kii'r  +  icS'f + }  •1^3 

where  m  is  an  integer. 

Proof:  The  proof  is  straightforward.  We  note  that  all  terms  are  square  integrable 
because  we  assumed  that  the  strain  energy  is  finite.  □ 

Let  us  also  note  that  if  the  strain  energy  is  finite  then  u  is  not  necessarily  square 
integrable  on  x  {—d,d)  not  only  for  the  possible  presence  of  rigid  body  motions. 
Nevertheless  all  terms  in  the  expressions  for  the  strain  energy  are  integrable. 

Theorem  3.1  allows  us  to  find  the  solution  of  the  3-dimensional  formulation  as  well 
as  that  of  any  n-model  explicitly  in  Fourier-transformed  form. 

Using  the  Fourier  transformation  of  generalized  functions,  we  now  write  down  so¬ 
lution  formulas  for  the  beam  and  plate  problems  on  the  unbounded  domain  floo.  We 
start  with  the  beam  problem  where  ^  ^  R}  because,  as  will  be  seen  in  Theorem  3.3, 
the  beam  solution  yields  immediately  the  plate  solution.  It  is  not  hard  to  prove: 

Theorem  3.2 


1.  Let  q(xi)  €  L2{R})  have  compact  support  and  assume  that 


/OO  eOO 

q{x)dx  =  /  xq{x)dx  =  0. 

•OO  */  — OO 


Then  the  Fourier  transform  of  the  beam  solution  —  (uj,U3)  introduced  in  Sec¬ 
tion  2  is 


where 


'^‘=1  Ik:":!  ]=[&’]’(«)■ 

^(^-^3)  =  {sinh(^3r3)  [(*/ -  i)  f  +  f  cosh(^X3)a}  (3.14) 

x(^>3:3)  =  ^  {cosh(^Z3)  [(1  -  i,)  2  +  |6]  -  ^sinh(^X3)a} 


a{(,d)  = 


cosh(^(i) 

cosh(^(i)  sinh(^d)  —  (d 


bU,d)  = 


sinh(^d) 

cosh(^(i)  sinh(^(/)  —  ^d 


and 
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The  solution  is  unique  modulo  the  space  of  rigid  body  motions  which,  under 
Fourier  transformation,  becomes  with  the  constraint  (2.8) 

=  {^*(6,2:3)  =  («i(6,a;3),-U3(6,a:3))  i  uj(6,2:3)  =  aa:3^(6), 

p\C 

“3(61^3)  =  ^^(6),  a,be  R} 


where  ^(^i)  is  the  Dirac  distribution  (in  ^i). 

2.  Let  q{x]  =  <fo  G  ,  ^0  ^  0.  Then  the  2r^Q^ -periodic  solution  of  the  beam 

problem  is 

Ui(xi,X3)  1  _  [  (oH^o,X3) 


u»  = 


0*2^1  ^0 


U3(^U^3)  J  [X(^0,l3) 

and  it  is  unique  up  to  the  rigid  body  motions 

Aper  =  {(ui,«3)  I  «i  =  0,U3  =  a,a  6  R^}. 


(3.15) 


Remark  3.1  In  the  special  case  of  periodic  solutions  considered  in  Theorem  3.1  we 
did  not  write  the  solution  u*’  in  its  Fourier  transformed  form.  Using  the  Fourier  trans¬ 
formation,  (3.15)  can  also  be  written  in  the  form  (3.13)  with 

9(0=2jr5(e-h6). 

Remark  3.2  The  distributions  ui,  U3  are  in  general  not  integrable  since  they  have 
a  singularity  at  the  origin.  Nevertheless  it  is  easy  to  see  that  their  strain  energy  in 
Theorem  3.1  is  finite  and  hence  the  solution  expressed  by  (3.13)  is  the  Fourier  transform 
of  the  unique  (up  to  rigid  body  motions)  solution  introduced  in  Section  2.  For  more 
details  see  [45]. 

Remark  3.3  In  Theorem  3.2  we  have  assumed  a  very  special  periodic  load  q.  It  is 
clear  that  we  could  represent  a  general  periodic  load  by  superposition. 

So  far  we  have  dealt  only  with  the  beam  problem.  Now  we  will  show  that  the 
Fourier  solution  formulas  for  the  beam  and  plate  problems  are  essentially  the  same. 
More  precisely,  we  have 

Theorem  3.3 


1.  Let  q{xi,X2)  €  L^iR^)  have  compact  support  and  satisfy  (2.9).  Then  the  Fourier 
transformation  of  the  solution  of  the  plate  problem  is 


u'’(6,6,a;3)  =  V5jr«7(6>^2) 


6  0(1^1. 2:3) 

60(1^1,2:3) 

X(I6,2:3) 


(3.16) 


M 


where  4^,  x  o,re  as  in  Theorem  3.2  and  the  space  of  rigid  body  motions  is  now 
given  by 

Afp  =  {u  =  (Ui(6,6,^3),U2(6,6,^3),U3(6,6-a;3))  I 

Ul  =  01X36(^1, 42),t^2  =  C2l3^(6,6),“3  =  6(5(^1,6)  -  Cii-7^  -  C2i-r^ 

0^2  C/si 


2.  For  q{x)  =  ^  ^  ^  0  we  have 


’  “5(6, 6, 2:3)  ‘ 

=  v/2ff 

rfi.oi/’(({;,o +&)■'", 5^3))  1 

= 

“5(6, 6, 2:3) 

.  “5(6, 6, 2:3) . 

Lx((^^.o+ao)‘'^3^3))  J 

and 


0+1242, 0) 


(3.17) 


=  {^^(6,612:3)  =  (ui,ti2,«3)  1  ui  =  0,U2  =  0,  U3  =  a,a  6  R^}. 


Proof:  First  we  observe  that  the  functions  t/>  and  x  given  in  (3.14)  are  even  in  ^  and 
hence  the  expressions  t/)(j^j,X3),  x(l^l»®3)  have  a  natural  meaning.  Let  us  now  consider 
the  expression  for  the  strain  energy  of  the  plate  given  in  (3.9).  Upon  the  substitution 


Ul 

U2 


(6, 6, 3:3)  = 


'6  6 

“1 

.  ^2  -6  . 

“2 

(656,2:3) 


^3(6, 6,  2:3)  =  \/^V3(6,6,2r3) 


we  get 


WITT)  +  +  + 


1(1  +  //) 

+  Y~^  M6wi  +  V3I*  d6d6dx3 

+(27r)-i 


L , .  ,(i^n^2i'+i^2n‘^6^6dx3. 


2{l  +  I/)  JB?x{-d4) 


Comparing  this  with  (3.10)  and  since  Q^i^)  is  given  by  (2.6),  we  get  U2  =  0  and  obtain 
from  Theorem  3.2 


Vl 

V3 


(656,2:3)  =  9(6,6) 


16^(16,2:3) 

x(I6,2:3) 


15 


Hence  we  have 


Ux 

6 

6 

U3 

"  lel 

.^2 

-^1  . 

V3 

and  this  implies  the  statement  of  the  theorem. 

In  the  case  of  the  periodic  problem,  i.e.  for  a  load  of  the  form  g  =  ^^e 

utilize  the  fact  that  q  is  the  Dirac  distribution.  □ 

In  what  follows  we  will  assume  that  and  are  given  by  (3.13),  (3.14),  (3.15), 
(3.16),  (3.17),  i.e.  we  are  fixing  the  rigid  body  motion. 


3.2.  The  Fourier  Transformation  of  the  Solution  of  n-models 


As  we  pointed  out  in  Section  2,  our  main  interest  will  focus  on  the  lowest  models  of 
the  hierarchy.  Because  the  same  relation  between  the  Fourier  transformation  of  the 
beam  and  plate  problem  holds  also  for  the  n-models,  we  will  only  address  in  detail  the 
expression  for  the  beam  model. 

For  the  (1,1,0)  model  we  use  B  =  z=  B{\\,  fix,  n)  where 


5(Ai,/xi,k) 


Ai  -f-  2/ii 

Aj 

Ai 

0 

0 

0 


Ai 

Aj  -f-  2//j 

Ai 

0 

0 

0 


Ai 

Aj 

A,  +  2/ii 
0 
0 
0 


0  0  0 

0  0  0 

0  0  0 

/ii  0  0 

0  K/ii  0 

0  0  K/Xi 


with  (see  (2.2),  (2.4)) 

. _ Eiux  _  Ex 

(1  +  i^i)(l  —  2i'i)’  2{l  +  ux) 

_  i'  p  _  E{1  +  2v) 

1  +  i/’  (l  +  t/)2 


(3.1S) 


(3.19) 

(3.20) 


and  K  >  0  is  arbitrary. 

For  the  n-models  n  =  {nx^n^iTij)  with  n  >  1,  nj  >  1,  03  >  2  we  use  B^  =  B(A,  //,  k) 
where  A,  //  are  given  by  (2.4)  and  k  >  0  is  again  a  free  parameter,  the  so-called  shear 
correction  factor.  Now  we  have  the  following  representations  of  the  solution. 


Theorem  3.4  For  every  n  and  B^  as  defined  above,  there  exist  functions 


SO  that  the  solution  of  the  n-model  admits  the  representation  (3. 16), (3. 17).  In 
particular,  we  have 
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1.  for  the  (1,0)  model 


and 


3i(l  —  v) 

(3.22) 

3(1 -t/)  2  1 

<Pi*  J  ’ 

(3.23) 

2.  for  the  (1,2)  model 

^  4Grfi  (30  +  rf!ic(l -i-)|fP)|{|'‘  ’’  '  ' 

3  30/c(l  —  I/)  +  €P1^|^(20  +  —  5ki/  —  2k'^i/)  —  d*K{l  —  i^)li^|'* 


AG(P  K(30  +  rf2^(l  _  j,)|^|2)|^|4 

2  15  3t/  +  <f^(i/-  l)|^p 

+X^ - 


3.  for  the  (3,2)  model,  as  |^|  —*  0, 


(3.25) 


3 

i;<p\  e 


(2  +  3kv)(P\ 


) 


X3  - 


2  —  KV 


6ki 


-x^  +  0(rf)  (3.26) 


+ 


(8  —  3kv)<P 

10k^2 


and,  as  |^|  — ►  oo, 


<f‘*(l  —  i/)(21k^(1  —  2i/)  +  (2  —  Ki/)^) 
525k2(1  _2i/) 


+  0(«f^e^),(3.27) 


(3.28) 


Proof:  The  formulas  (3.22)-(3.28)  are  derived  exactly  in  the  same  way  ais  the  formulae 
(.3.14),  (3.15).  □ 

Let  us  remark  at  this  point  that  the  functions  ^  and  x  *^re  meromorphic  in  ^  in  a 
neighborhood  of  ^  =  0. 

The  main  idea  is  now  to  analyze  the  error  of  the  n-model  in  the  Fourier  transform 
variables.  We  note  that  the  matrix  5“  is  defined  differently  for  723  =  0  and  >  2. 
In  both  cases  we  included  the  shear  factor  k  and  will  analyze  its  optimal  value  which 
leads  to  the  minimal  error  when  the  solution  of  the  n-model  is  compared  with  the  three- 
dimensional  solution.  We  could  also  determine  all  the  coefficients  in  by  the  same 
principle.  This,  however,  will  lead,  for  the  materials  considered  here,  to  our  matrices 
B  with  K  being  the  only  free  parameter. 
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Remark  3.4  The  introduction  of  the  shear  factor  is  a  very  old  idea.  It  was  introduced 
for  the  dynamic  plate  problem  as  well  as  for  the  static  one.  In  [2,  40,  46]  the  shear 
factor  K  =  5/6  is  recommended.  In  [21,  32,  51]  a  value  for  the  shear  factor  is  obtained 
so  that  the  velocity  of  short  surface  waves  is  modelled  accurately.  In  [13,  16]  the  shear 
factor  is  designed  so  that  an  accurate  solution  will  be  obtained  for  a  cantilever  beam 
with  a  given  cross  section.  Various  other  shear  factors  were  suggested  for  laminated 
plates,  e.g.  [10,  35,  49,  50].  In  [36],  an  adaptive  construction  for  the  shear  factor  is 
proposed.  All  these  suggestions  are  ad  hoc.  For  the  (1, 1,0)  model  [28,  29,  38]  analyze 
the  optimal  selection  of  the  shear  factor  in  connection  with  nonhomogenous  lateral 
boundary  conditions  and  <7  =  0  (the  so-called  Saint- Venant  problem  for  plates)  and 
show  that  for  a  particular  value  of  k  the  energy  error  is  asymptotically  of  higher  order. 

In  contrast  to  these  derivations  we  wilt  present  in  the  next  section  an  approach 
which  leads  to  uniquely  determined  special  values  of  the  shear  factor  which  in  turn 
give  a  higher  asymptotic  rate  of  convergence  for  the  data  of  intent  in  the  plate  bending 
problem. 

4.  Error  Analysis  and  Optimal  Shear  Factors  for  the  Hierar¬ 
chical  Models 

The  Fourier  transformation  allows  us  to  analyze  the  performance  of  the  n-models  with 
respect  to  various  measures  of  interest.  In  this  section  we  will  address  some  of  them. 

4.1.  The  Load  Spaces  and  the  Basic  Performance  Measures  “• 

In  Section  2.2  we  have  introduced  the  accuracy  measure  relatively  to  a  set  of  loads  V. 
Let  us  now  define  a  specific  set  V  =  j  of  loads  g{x\,X2)  for  the  plate  problem  and 
g{xi)  for  the  beam  problem.  Let  for  <  >  1,  f  >  0  integers  and  s  €  i?’, 

B',.,  =  {?u(i)  =  o,  for|i|>/c,/(i  +  iep)’i?(e)i"<i{<i, 

(D“9)(0)  =  0,  q  =  (01,02),  0  <  0i,02  integers  01+02  =  |o|,  |o|  <  <, 

/  =  A,-  >  0,  j  =  0, 1, . . . ,  /}. 

We  note  that  from  the  bcisic  properties  of  the  Fourier  transformation  of  distributions 
with  compact  support  we  have  for  |^|  <  1,  jo]  <  s  that 

KO'«)«)l  <  c(c.,»)  (/  |jp(i  +  |«p)-<;«)'''.  (4.1) 

and  hence,  because  t  >  1,  we  obtain 

/  <  00. 

4|<l<> 


The  space  Sj  ,  is  a  space  of  generalized  functions  (distributions),  based  on  which  we 
introduce  now  the  basic  performance  measures  E  corresponding  to  the  data  of  interest 
Q.  Obviously  we  could  introduce  many  others.  We  have 


a)  The  energy  error  performance  measure. 

b)  The  midsurface  deflection  measure. 


Emid{B,  n,  q) 


|/«.lx(|{|.0) -V(i;i.0)l=l« 


(4.3) 


c) 


The  average 
Let 

and 


deflection  measure. 

1 

U3.AV(0  =  ^U3{(,X3)dx3 

^”av(0  =  ^3)^x3. 


Then 

EAv{B,n,q)  = 

The  following  theorem  shows  that  in  order  to  analyze  the  energy  performance  measure 
asymptotically  as  d  — >  0  we  need  only  consider  the  numerator  of  (4.2). 

Theorem  4.1  Let  q  €  B[  „  s  >  —1/2,  /  >  0,  0  <  d  <  do  <  1-  Then 


IfelxAvdfD-^-SvdeDPljP^'^ 
l/„.  LxAvdfl)Pl9l'<i«r'' 


/p  .xdfU)l«OPrf{  >  (4.5) 

where  C  >  0  is  independent  of  d. 

Proof:  First  we  observe  that  it  follows  from  (3.14)  that  for  \^\d  >  1 
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with  Cl  independent  of  d.  Further  we  have  for  all  \^\d  <  1 

x{\^,d)  =  aod-^\(\-^-\-RU,d)  where  d)l  <  (4.7) 

and  Cj  is  independent  of  d  and  oq  ^  0.  Now  we  write 

I  Xm,d)\q\^di=  f  +/  +/ 

jR^  Jm>i  yKi<i  J\i\>u\i 


and  consider  the  three  terms  separately.  Using  (4.6)  we  get  for  s  >  —1 12  that 

/  x(kU)i9p^^  <  I/.  x(ki,^)(i  +  i^i^r(i  +  i^i'ri9i'‘^^ 


<  sup  |x(|{|,<i)(i  +  ieh-i<c. 


Further 


4,  x(\Mm(  =  -0^-^  4.  + X,,,  m 


'l€l<i 

As  q  €  Bi  j,  we  have 

and 

Hence  for  d  <  do  we  get 
Finally 


aod-^  f  \(\-^\q\^d^  =  aod-^Ao 

1/  R{^,d)\qM 

\J  t\<l 


<  d-^C. 


f  =  a^d  ^(/lo  +  0(d^)). 
dKI<i 


where 


/  xm,d)\q\^d^  =  aod-^  f  \^\-^\ed^+l  R{i,d)\q\^d^ 

aod~^  [  =  aod~^p,  p>0 

J\i\>iM\d<l 


and 


/  R{Cd)\qM 


1 


<  cd-’  j  . ^(1  +  i^n-*(i  + 


i^p 


<  Cd-\ 
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Hence  collecting  all  the  terms  we  get  for  0  <  <  £?o 


>  |ao|.i-=(l  +  0(<i")) 


which  yields  (4.5).  □ 

In  what  follows  we  will  also  need  the  behaviour  of  XiXav  for  large  |^|.  We  have 
Theorem  4.2 


For  l^ld  >  1, 

iX>lv(l^|)|  <  |^|2^ 

(4.8) 

lx(l^l,0)l  <  m>0  arbitrary. 

(4.9) 

For  |^|d  <  1, 

X^v(^,d)  =  aod-^|^r''  +  7?(f.d) 

(4.10) 

where 

|/^l<Cd-M^|-^ao^0, 

and 

x(iao)  =  5od-^iei-‘‘  +  /?, 

i/2|<Cd-^lei-',do7^0. 

(4.11) 

Proof:  (4.9)  and  (4.11)  follow  directly  from  (3.14). 

By  an  easy  calculation  we  obtain  from  (3.14)  also  that 

XAviO  = 

1  (3  -2i/)(l  -  +  4^de2«'^ 

G  4^2</(i  _e-««<i)  + 

from  which  (4.10)  and  (4.11)  follow. 

Now  we  can  prove  along  the  lines  of  Theorem  4.1  the  following  result  for  \av  and 

xK.o). 

Theorem  4.3  Let  g  €  B3 ,,  /  >  0,  s  >  —2  and  0  <  d  <  do  <  1.  Then 

J  \xAy{0?\q?df>Cd-\  (4.1.3) 

and  for  every  real  s  we  have 

j  ix(«.o)riiiV{ 

where  C  >  0  is  independent  of  d. 


(4.12) 

□ 


(4.14) 
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Theorems  4.1  and  4.3  allow  us  to  analyze  only  the  numerators  in  (4.2),  (4.3)  and 
(4.4)  in  order  to  prove  the  conditions  (Hi)  and  (H2)  in  Definition  2.2.  Since  the  analysis 
of  all  error  measures  proceeds  along  similar  lines,  we  present  it  here  only  for  (4.2)  in 
detail  and  elaborate  briefly  on  (4.3),  (4.4).  The  following  theorem  shows  that  the 
asymptotic  behaviour  of  the  modelling  error  for  the  n-model  is  completely  governed  by 
the  behaviour  of  near  |^|  =  0,  provided  the  data  is  sufficiently  smooth. 

Theorem  4.4  Assume  that  for  |^jd  <  1  and  every  integer  k  >  0 


xm,d)  = 

i=o 


(4.15) 


and  analogously  that 


<  C-d-3+2(fc+l)l^|-4+2(fc+l) 


(4.16) 


V(kl,d)  = 

j=0 


(4.17) 


<  (7</-3+2(t'+l)|^j-4+2(/:+l)_ 

In  addition  assume  that  for  |i^|d  >  1 


(4.18) 


<<?;- 


(4.19) 


where  C  is  independent  of  d.  Let  further 

r  =  mjn{aj  ^^a^] 

where  the  minimum  is  taken  over  all  nonnegative  integers. 

Then  for  q  €  B[ ,,  s  >  —1/2,  I  >  j'  we  get  that 

Vim  EE{B,n,q)d~^  =  C,  0  <  C  <  oo  for  /3  =  min{s  +  2,i*}. 

d—*0 

proof:  Using  Theorem  4.1  we  have  to  prove  that 

=  Cd^^-^(l+oil)) 


(4.20) 


(4.21) 


as  d  —y  0.  To  this  end  we  w.ite  as  in  the  proof  of  Theorem  4.1 


z  =  I  +1  +/ 

Now  using  (4.19)  we  get  for  s  >  —1/2 

[  (,v(l«U)-"'x"(ifi.<i))i9p<;{  <  sup  c(i  +  ifr)-i 


l£l.l£l<i>i 


Further, 


^l£l<i.l£M<i 
where  \R\  <  and  finally 


f  =  (aj.  -  V“  )d-3+2j‘  Aj.  +  R 


f  =  {aj.  +  R 


□ 


where  |/i|  <  Cd  3+20*+i)  and  therefore  (4.21)  is  proven. 

Analogously  we  can  prove 

Theorem  4.5  Suppose  that  for  XAvi\^\),  x(l^l»0)  and  ®x”(kl>0)  have 

expansions  analogous  to  (4.15)-(4.18)  for  <  1  and  that  for  |^|d  >  1  there  holds 


h..u(l41)-®x5,.(l«l)| 

|x(l«l,o)  -®x"(l«l.o) 


< 


< 


c 

m' 

c 

m' 


(4.22) 

(4.23) 


Then,  if  q  E  with  s  >  —2  and  I  >  2j"  —  2,  we  have 


\\m  EAv{B,n,q)d  ^  =  C  0  <  C  <  oo,  with  /3  =  min{s  +  4, 2/*}.  (4.24) 

d— *0 


Here  j‘  is  defined  as  in  (4.20).  For  the  error  measure  E\riD  defined  in  (4.3)  the  same 
bound  as  in  (4.24)  holds. 

The  proof  Theorem  4.5  is  analogous  to  that  of  Theorem  4.4  and  uses  Theorem  4.2. 

From  Theorems  4.4  and  4.5  we  see  in  particular  that  the  n-model  is  asymptotically 
exact  only  if  j*  >  1. 
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4.2.  Asymptotic  Analysis  for  n-Models 

We  shall  now  apply  theorems  4.3  and  4.4  to  the  analysis  of  the  n-rrodels  and,  in 
particular,  to  the  determination  of  the  optimal  shear  correction  factors  k. 


4.2.1.  The  Reissner  -  Mindlin  —  (1,1,0)  Model 
We  have  from  Theorem  3.2  that  for  |^|d  <  1 


X(ki,a:3) 


1  /3(1-//)  Zuxl  3(8 -3j/) 

AG  V  10|^|2d 

(157i/  -  227)d'‘  +  (1050  -  630i/)d^x2  _  175(1  ^ 

1400d3 


(4.25) 

+  0{d^\^f) 


and 


Xav(I^I)  = 


J_ 

AG 


/3(l-i/)  12 -7i/  lld(l-i/)\ 

I  175  ) 


+  o{d^m. 


(4.26) 


On  the  other  hand  using  Theorem  3.4,  we  get  for  n  =  (1, 1,0) 


3(1 -t^)  2  1 


(4.27) 


By  (4.7)  and  (4.27)  we  see  that  (4.19)  is  satisfied  and  furthermore,  comparing  (4.27) 
with  and  (4.25),  we  find  for  the  mecisure  Eg  that: 

If  (1  G  3  for  /  >  1  and  s  >  —1/2  then 


0  = 


min{2,  s  +  2}  if  «  =  ^opt  = 

1  if  K  ^  ^opl 


For  the  measure  £mid  we  obtain  with  Theorem  4.5  in  the  same  fashion: 
If  7  G  ^3  3  for  /  >  2  and  s  >  — 2  then 


(4.2S) 


f  min{4,5  +  4}  if  ;c  =  /Cop^  =  3,^ 

(2  if  «  7^  /Copt 

For  the  measure  jFav  we  get: 

If  <7  G  B3  3,  /  >  2  and  s  >  —2,  then 


min{4,s  +  4}  if /c  = /Copt  = 
2  if  /C  ^  Hapt 


(4.29) 


(4.30) 


Let  us  comment  on  these  results.  We  see  that  the  optimal  value  of  k  depends  on 
the  goal  of  the  plate  modelling.  For  =  0  we  get  in  all  cases  that  «;„pt  =  5/6  i.e. 
the  classical  value  proposed  for  example  in  [2,  40,  46].  Secondly,  for  data  q  which  is 
not  sufficiently  smooth  (or  equivalently,  for  the  case  when  the  solution  of  the  three- 
dimensional  problem  lacks  regularity),  i.e.  when  s  is  small,  nothing  can  be  gained  by 
the  above  choices  of  /c.  We  will  analyze  this  case  seperately  below.  Finally  we  remark 
that  for  the  choice  B  =  A  the  (1, 1,0)  model  is  not  asymptotically  exact,  i.e.  the  use  of 
the  modified  material  matrix  in  (3.18)  is  necessary  for  the  asympf-^tic  exactness 

of  the  (1,1,0)  model. 


4.2.2.  The  (1, 1,2)  Model 

From  (3.26),  (3.27)  and  (4.26)  we  find  estim.ates  which  are  identical  to  (4.2S)-(4.30) 
provided  that  the  optimal  shear  correction  factor  is  given  by 


12 -2u 


^opt  - 


-1  + 


1  + 


201/2 


(12-2i/)2 


1/  ^  0. 


(4.31) 


We  point  out  that  unlike  for  the  (1,  l,0)-model  here  /Copt  is  the  same  for  all  three 
accuracy  measures  used.  In  Table  1  we  report  the  optimal  shear  correction  factor  for 
various  values  of  v. 


Table  1  The  values  of  /Copt  for  the  (1,1,2)  model. 


V 

^opt 

0 

5/6 

0.25 

0.S67 

0.30 

0.S74 

0.40 

0.S87 

0.50 

0.899 

For  other  values  of  the  Poisson  ratio  i/  it  could  be  easily  computed  from  (4.31). 

Remark  4.1  For  /c  =  1  the  measure  Ee  is  also  the  relative  error  in  the  energy  norm. 
For  K  ^  I  the  error  in  the  energy  norm  has  to  be  computed  separately  and  we  shall 
not  elaborate  on  it  here. 


4.2.3.  The  (3,3,2)  Model 

Reasoning  as  for  the  (1,1,2)  model  we  obtain  now  that 
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For  the  error  measure  Ee'-  If  9  G  with  s  >  —\l'2  and  I  >  2,  then 

n /  min{2,  ^  2}  if  /c  —  /^opt  —  1) 

^  ~  I  1  ii  K  ^  «opt 

For  the  error  measures  Eav  and  Em  id  we  find:  If  ^  G  with  s  >  —2  and  /  >  2,  then 

o_f  min{4,5  +  4}  if  #c  =  Kopt  =  1, 

^  ~  [  2  if  /£  7^  /Copt 

Consequently,  for  sufficiently  smooth  loads  q,  the  introduction  of  shear  correction 
factors  K  1  into  the  n-models  with  n  >  (3,3,2)  componentwise  will  always  decrease 
the  rate  of  convergence  and  hence  the  choice  =  /I  is  optimal  in  these  cases.  We 
further  see  that  the  rate  of  convergence  /3  achieved  by  the  (3,3,2)  model  is  the  same 
as  for  the  (1,1,0)  model  with  /Copt-  For  sufficiently  large  s  we  get  a  higher  rate  of 
convergence  only  by  using  the  (3,3,4)  model.  Then  we  can  replace  2  and  4  in  the 
expressions  for  /?  by  3  and  6,  respectively  (see  [45]). 


4.3.  The  Error  Analysis  for  Unsmooth  Loads 

In  the  previous  sections  we  analyzed  the  behaviour  of  the  plate  modelling  error  as 
d  — >  0  for  a  class  of  loads  q.  VVe  have  seen  that  the  accuracy  of  the  n  model  was  in  this 
case  governed  by  the  behaviour  of  the  functions  0  and  \  for  small  <f. 

Furthermore,  the  smoothness  of  the  load  q  in  the  class  B[^  was  characterized  by  the 
parameter  s  which  in  turn  characterized  the  decay  of  \q\  as  |if)  — >  00.  In  the  proof  of 
Theorem  4.3  and  Theorem  4.5  the  minimal  allowable  s  was  governed  by  the  decay  of 
(  y  — for  large  Therefore  we  could  be  interested  in  the  selection  of  the  shear 
correction  factor,  for  which  this  decay  will  be  maximal. 


Theorem  4.6  Consider  the  (1,1,0)  model.  Then  we  have  for  |^]d  >  1  that 


|xav(0  -®Xav(0|  <  77^  forK  =  - 


3-2i/ 


=  K, 


opt 


Cl 


^  |xAv(‘f)  ~^Xav(0|  ^  lf|2J  ^  ^ 


(4.33) 

(4.34) 


^j^j2  |A«.  V-,/  AAVV>/|  -  |^J2J 

where  C,  C\,  C2  are  positive  constants  which  are  independent  of  d. 

« 

Proof:  Using  (4.12)  and  (4.27),  it  is  easy  to  show  that  there  exist  Ci,C2  >  0  so  that 


c, 


1 


< 


XAviO  -®Xav(0  ~  9  ( 


i  -1]  — 


<C2- 


1 


kJ ed\-  \f\^cp 


from  which  (4.33)  and  (4.34)  follow.  □ 

This  theorem  can  now  be  used  to  obtain  optimal  shear  correction  factors  for  the 
error  measure  Eav  and  certain  classes  of  data  q.  To  this  end  let  us  assume  that  the 
thickness  d  is  fixed.  Then,  based  on  Theorem  4.6  we  can  construct  a  sequence  of  leads 
Qk,  k  =  1,2, ...,  so  that 

\XAv\'^\qk?d^  <  oo, 

I  oo  as  k oo 

J 

with  qk  €  where  I  >  0  and  — 4  <  s  <  — 2  and  for  which 


2 

EAviB,n,qk) Q  for  k  =  Kopt  =  ■— 


and 

EAv{B,n,qk)  0  for  k  ^  Kopt- 

In  Table  2  we  show  the  values  Wopt  from  (4.30)  for  “smooth”  q  and  (4.33)  for  “un¬ 
smooth”  q  for  the  (1, 1,0)  -  model  and  the  error  measure  Eav- 


Table  2  Comparison  of  «opt  of  (4.30)  and  (4.33)  for  Eav  and  the  (1,  l,0)-model. 


u 

(4.30) 

(4.33) 

0 

0.83 

0.67 

0.25 

0.98 

0.80 

0.30 

1.01 

0.83 

0.40 

^T09 

0.91 

0.50 

1.18 

1.00 

Table  2  shows  that  in  the  case  of  average  displacement  we  take  for  less  smooth  q  the 
smaller  shear  factor.  We  emphasize  that  this  conclusion  is  valid  only  for  the  average 
normal  deflection  measure  E^v- 

Remark  4.2  We  point  out  that  for  the  measures  Emid  and  Ee  the  optimal  shear 
correction  factors  Kk  for  the  above  sequence  of  loads  must  become  large  and  positive 
cis  «:  — >  oo.  This  follows  from  a  comparison  of  (4.27)  and  (3.14). 

In  section  3  we  have  introduced  the  periodic  solution  with  the  load  q  =  e'^^,  i  = 
(xijZj).  In  this  case  the  solutions  of  the  three  dimensional  problem  and  the  n-model 
have  the  same  form  in  (ii,X2).  Consequently,  we  can  eaisily  compute  the  shear  factor 
K  as  a  function  of  |^o|d  which  yields  the  same  result  for  the  exact  and  n-model  solution 
for  the  given  quantity  of  interest,  for  example  for  the  measure  Eav-  The  graphs  of  the 


Figure  1:  Optimal  shear  correction  factors  for  the  (1,1,0)  model  (top)  and  the  (1,1,2) 
model  (bottom)  with  periodic  loads  as  functions  of  |^old. 
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optimal  shear  correction  factors  for  the  (1, 1,0)  and  (1, 1,2)  model  are  given  in  Figure 
1.  We  note  that  in  the  limit  |^o|d  0  we  obtain  from  (4.30),  while  for  |^o|d  oo 

we  get  Kopt  from  (4.33). 

In  summary,  we  saw  in  this  section  that  the  shear  correction  factor  for  any  n-model 
can  be  uniquely  determined  by  some  optimality  requirement  and  we  have  analyzed  some 
of  the  possible  choices  for  this  requirement.  We  found  that  above  a  certain  threshold  of 
the  model  order,  the  introduction  of  a  shear  correction  factor  into  the  plate  model  docs 
not  improve  the  asymptotic  rate  of  convergence.  For  the  homogeneous  and  isotropic 
materials  considered  here  only  the  (1, 1,0)  and  the  (1, 1,2)  models  can  be  improved  by 
a  judicious  choice  of  /c.  In  addition,  the  optimal  shear  correction  factors  also  depend 
on  the  regularity  of  the  three-dimensional  solution  which  is  approximated  by  the  plate 
model.  We  point  out  that  many  error  measures  different  from  the  one  considered 
here  could  be  analyzed  in  the  same  fashion.  It  is  interesting  to  note  that  for  different 
goals  of  the  plate  modelling  as  well  as  for  different  regularity  properties  of  the  three 
dimensional  solution  different  shear  factors  are  recommendable.  This  also  suggests  that 
in  the  neighborhood  of  the  boundary  different  shear  factors  should  be  used  than  in  the 
interior  since  boundary  layers  or  singularities  of  the  three  dimensional  solution  govern 
the  quality  of  the  plate  model  there.  We  will  also  see  this  in  the  numerical  experiments 
in  the  following  section. 

5.  Numerical  Experiments 

In  the  previous  sections  we  obtained  several  optimal  optimal  values  for  the  shear  cor¬ 
rection  factor  K  for  the  infinite  beam  and  plate  problem  and  for  different  quantities  Q 
of  interest. 

For  the  problems  in  applications  the  domain  is  bounded.  In  this  case  the  theoretical 
convergence  of  the  various  models  can  be  analyzed  as  in  [7]  where  the  model  (1, 1,0) 
was  addressed.  Then  also  the  solution  is  usually  not  smooth  because,  for  example, 
of  various  boundary  layers  caused  by  the  incompatibility  of  the  data  and  the  lateral 
boundary  conditions,  see  e.g.  [4,  5].  Hence  we  will  now  investigate  the  influence  of  k. 
for  the  following  model  problem. 

Consider  a  clamped  beam  of  thickness  2d  =  0.2  and  length  2  with  Youngs  modulus 
E  =  I  and  Poisson  ratio  u  =  0.3  as  in  Figure  2  where  the  load  q{i)  is  also  defined. 
Because  of  the  obvious  symmetry  in  the  problem  we  do  the  computation  for  the  half 
beam  only.  First,  as  exact  solution  we  will  use  the  numerical  solution  by  the  h-p  finite 
clement  method  obtained  with  the  program  MSC/PROBE  with  high  polynomial  degree 
p  =  S  and  refined  meshes  near  Xi  =  0.2,  Xi  =  0.5  and  13  =  ±d.  The  relative  error  in 
the  energy  norm  was  less  than  10~®  so  that  this  solution  could  be  considered  as  the 
exact  one. 
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Figure  2:  Scheme  of  the  model  problem 

Then  the  n-model  of  the  beam  problem  leads  to  a  two  point  boundary  value  problem 
for  a  system  of  ordinary  differential  equations.  This  boundary  value  problem  was  in 
turn  solved  numerically  to  high  accuracy  by  the  method  described  in  [3]  and  also  its 
solution  could  be  taken  as  the  exact  solution.  Hence  the  difference  between  the  two 
solutions  is  the  modelling  error  as  discussed  in  the  previous  sections. 

Figure  3  shows  the  pointwise  error  in  the  average  normal  deflection  for  the  (1,0)- 
model.  The  curves  correspond  to  (from  above)  «  =  1.0,  0.938  (the  theoretically  optimal 
value  for  the  midsurface  deflection  from  (4.29)  0.935,  0.87  (suggested  in[13])  and  k  = 
5/6  which  was  suggested  In  [2,  40,  46|. 


Figure  3:  Pointwise  error  in  the  average  normal  deflection  for  the  (l,0)-model. 

(1)  K=l 

(2)  K  =  0.938  (optimal  value  (1.29)  for  the  inidfiber) 
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(3)  K  =  0.935  (optimal  numerical  value) 

(4)  K  =  0.87  (Cowper  [13]) 

(5)  K  =  5/6  (Reissner  [40]) 

Clearly  the  results  obtained  with  «  =  1  exhibit  a  large  error  which  could  be  expected. 
Surprisingly,  however,  the  values  proposed  in  [13]  and  [2,  40,  46]  are  not  leading  to 
a  smaller  error,  either.  On  the  other  hand,  the  value  k  =  0.935  which  was  found 
numerically  as  an  optimal  value  for  the  given  case  is  very  close  to  the  value  k  =  0.93S 
which  is.  according  to  (4.29),  optimal  for  the  midsurface  deflection.  VVe  point  out  that 
K  =  0.935  is  not  close  to  the  value  k  =  1.01  in  (4.30)  respectively  k  =  0.83  in  (4.33). 
The  optimal  value  is  obviously  in  between.  This  effect  is  caused  by  the  fact  that  the 
solution  in  not  smooth  uniformly  through  the  beam.  Consequently,  (4.33)  suggests 
the  value  k  =  0.833  near  xi  =  0.5  and  by  (4.30)  k  =  1.01  elsewhere.  Because  of  the 
boundary  layer  we  select  therefore 

/c(i)  =  1.01  +  (0.833  -  1.01)e-'°<^-“-')'.  (5.1) 

In  Figure  4  we  show  the  error  of  the  (l,0)-model  for  the  average  (through  the  thickness) 
normal  deflection  for  the  numerically  optimal  value  k  =  0.93  and  for  k{x)  given  by  (5.1). 


Figure  4:  Pointwise  error  in  the  average  normal  deflection  for  the  (l.O)-model. 

(1)  K  =  0.935  (optimal  numerical  value) 

(2)  K  =  k(x)  as  in  (5.1) 

We  see  that  this  x-dependent  shear  factor  performs  as  well  as  the  best  constant  value 
for  K. 
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So  far  we  addressed  the  (l,0)-model.  In  Figure  5  we  show  analogous  results  for 
the  (l,2)-model.  We  depict  the  error  in  the  average  normal  deflection  for  the  constant 
shear  correction  factors  «:  =  1.0,  k  =  0.874  from  (4.31)  and  for  the  variable  k 

k{x)  =  0.874  +  (0.833  -  0.874)e-^°('-°-5)'  (.5.2) 

with  the  value  k  —  0.833  for  the  nonsmooth  solution  (see  Figure  1).  We  see  once  more 
the  same  effect  as  for  the  (l,0)-model. 
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Figure  5:  Pointwise  error  in  the  average  normal  deflection  for  the  (l,2)-model. 

(1)  K  =  1 

(2)  «=  0.874  (optimal  value  (4.31)) 

(3)  K  =  k{x)  (variable  correction  factor  as  in  (5.2)) 

Remark  5.1  In  practice  we  would  select  the  shear  factor  piecewise  constant  over  the 
finite  elements  used  to  solve  the  beam/  plate  problem  approximately  rather  than  a 
function  which  depends  smoothly  on  x. 

Finally  we  show  in  Figure  6  that  the  value  of  Kopt  given  by  (4.31)  which  is  the  same 

for  all  criteria  used  is  also  in  practice  (i.e.  for  the  bounded  domain)  very  good  for 

various  accuracy  measures. 
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Fig  .e  6:  Error  in  various  quantities  in  dependence  on  k  for  the  (1,2)  model. 

(1)  4‘-='’(0,0)-U3(0,0) 

(2)  Error  in  the  average  normal  deflection  at  zi  =  0 

(3)  4''^^(0.0‘0-«3(0,0.4) 

(4)  Energy  of  the  (1,2)  model  -  energy  of  the  e.xact  solution 


6.  Conclusion 

We  have  presented  the  concept  of  hierarchical  modelling  as  a  numerical  method  for 
the  approximation  of  a  three  dimensional  elasticity  problem  on  a  “thin”  domain.  This 
method  can  be  ea.sily  implemented  in  the  context  of  the  h-p  version  of  the  finite  element 
method  (which  is  available  in  the  code  MSC/PROBE).  The  selection  of  the  shear 
correction  factors  Wcis  then  based  on  an  optimization  of  this  method.  The  approach 
that  we  have  presented  suggests  concrete  values  for  the  shear  correction  factors,  and  we 
found  that  the  optimal  choice  depends  on  the  aims  of  computation.  We  point  out  tha<^ 
the  optimal  shear  correction  factors  for  the  error  measures  considered  here  are  identical 
for  the  (1, 1,2)  model.  In  contrast,  for  the  Reissner  Mindlin  (1, 1,0)  model,  the  optimal 
shear  correction  factors  strongly  depend  on  the  error  measures  used.  We  further  showed 
that  the  optimal  shear  correction  factor  should  also  depend  on  the  smoothness  of  the 
solution  and  consequently  have  different  values  in  the  interior  and  in  a  vicinity  of  the 
boundary  of  the  plate. 

In  this  context  we  also  remark  that  the  Reissner  Mindlin  and  the  (1,1,2)  model  have 
quite  a  different  boundary  layer  behaviour  on  bounded  domains  [4], [5].  It  is  possible  to 
use  different  models  from  the  hierarchy  simultaneously  in  different  parts  of  the  plate. 
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for  example,  in  a  neighborhood  of  the  boundary  the  (1, 1,2)  or  a  higher  order  model, 
while  in  the  interior  the  (1,1,0)  model  with  the  proper  shear  correction  factor  could 
be  used. 

Although  we  presented  here  only  an  analysis  of  the  isotropic  case,  we  emphasize  that 
our  approach  can  also  be  utilized  for  anisotropic  and  laminated  materials  (where  the 
approximation  (2.20)  will  not  be  polynomial  in  I3  any  more)  as  well  els  for  the  analysis 
of  problems  other  than  the  classical  plate  bending  problem,  and  these  generalizations 
will  be  the  subject  of  a  forthcoming  paper. 
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